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Abstract. DNA renaturation is the recombination of two complementary single 
strands to form a double helix. It is experimentally known that renaturation proceeds 
through the formation of a double stranded nucleus of several base pairs (the rate 
limiting step) followed by a much faster zippering. We consider a lattice polymer model 
undergoing Rouse dynamics and focus on the nucleation of two diffusing strands. We 
study numerically the dependence of various nucleation rates on the strand lengths and 
on an additional local nucleation barrier. When the local barrier is sufficiently high, 
all renaturation rates considered scale with the length as predicted by Kramers' rate 
theory and are also in agreement with experiments: their scaling behavior is governed 
by exponents describing equilibrium properties of polymers. When the local barrier is 
lowered renaturation occurs in a regime of genuine non-equilibrium behavior and the 
scaling deviates from the rate theory prediction. 
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1. Introduction 

When a double-stranded DNA molecule is heated above some characteristic temperature 
it undergoes a denaturation or melting transition and the two strands separate (21]. 
The reverse process can also occur: two complementary single stranded nucleic acids at 
temperatures below their melting point bind to form a double-helix. This transition is 
known as renaturation [I]. 

The first investigation of the DNA renaturation kinetics dates back to the work 
of Wetmur and Davidson of the late sixties [22]. They experimentally measured the 
renaturation rate k 2 and found it scales as a function of L, the length of the strands, 

k 2 ~ L a (1) 

with a = 0.5. A recent detailed analysis of experimental data for DNAs from different 
organisms, with different conditions of salt and temperature was performed by Sikorav 
et al. [TTJ who estimated a = 0.51(1). 

The origin of the exponent a has been discussed in the literature [22l [T6l [T7] . 
Renaturation is believed to proceed via two steps: the formation of an active nucleus 
of a few base pairs, which is followed by a rapid zippering until the double helix is 
formed over its full length [I]. The rate-limiting step is the nucleus formation, therefore 
theoretical investigations of k 2 focus on the nucleation mechanism from the two separate 
strands. It was first believed that a = 1/2 could be explained by the theory of diffusion 
limited reactions [IB] . This is however at odds with more modern concepts of polymer 
physics, as explained in Ref. [17]: the correct analysis of diffusion limited reaction 
would imply a scaling of the type k 2 ~ L, which is not consistent with the experiments. 
According to Sikorav et al. [17], experimental data suggest that renaturation is an 
activated process. Using Kramers' theory, in which the polymer dynamics is mapped 
into that of an effective Brownian particle escaping a potential barrier, they find that 
the rate of nucleation of two monomers scales as [IT] : 

k mon ~ (2) 

where the exponent a 4 is that associated to the equilibrium contact probability through 
the middle monomer of the two renaturating polymers (see next section). As nucleation 
can occur in any of the L monomers of the two approaching strands, one expects [IT] 
the renaturation rate to scale as: 

k 2 ~ Lk mon ~ L 1+ °\ (3) 

Using the numerical estimate 04 = —0.48 [II] for self- avoiding walks in three dimensions 
one gets a = 1 — 0.48 = 0.52 in good agreement with experiments. 

The argument leading to Eq. ([2]) is based on a mapping of the DNA renaturation 
problem onto that of a single Brownian particle crossing a free energy barrier [IT] . The 
barrier has an entropic origin: it is generated by excluded volume interactions of the 
two approaching strands. One applies then Kramers' rate theory, which states that the 
transition rate is proportional to e AF / fcsT for a free-energy barrier AF and a temperature 
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T. For two polymers of length L one has AF = ksTa^logL (see Eq. (171)), from which 
one obtains Eq. (T5]). 

Mapping the complex polymer dynamics into that of a single effective Brownian 
particle performing a one dimensional motion is a huge simplification of the problem. 
Moreover Kramers' approach assumes that the polymer dynamics is sufficiently slow 
so that at any given time the conformation of the polymer is well approximated by its 
equilibrium distribution. 

To gain more insights into the origin of the exponent describing the scaling behavior 
of the rate constant ft2 as observed in experiments, we performed numerical simulations 
of a coarse-grained model of DNA renaturation. The model is similar to those used 
in the past for the study of equilibrium properties of the DNA denaturation transition 
[6j [5j |9]. It does not incorporate helical degrees of freedom of the double stranded 
DNA as the rate limiting step is nucleation, i.e. the first "reactive" contact between 
the approaching strands. The lattice nature of the model allows the simulation of long 
times and long polymers and, thus, an accurate determination of the scaling properties 
of the renaturation rates. 

The dynamics of the DNA denaturation has attracted quite some attention in the 
recent literature (3J [TUl [2] , but the reverse process of renaturation has been less studied. 
Our results show that the renaturation process has a rich dynamics and that the scaling 
predicted by Kramers' rate theory (Eqs. (T5j) and (J3J)) is only valid in a specific limiting 
case. 

This paper is organized as follows: We start in the next section with a review of a 
few theoretical results about scaling properties of equilibrium free energies. In Section [3] 
we introduce the model used in this study and in Section H] we present its equilibrium 
properties. The simulation results for the rate constants for renaturation are presented 
in Section [5j Finally, results are discussed in Section |6j 

2. Scaling properties of equilibrium free energies 

Before entering into the details of our simulations, we briefly recall the scaling properties 
of the free energy required to bring two self avoiding walks (SAWs) at a given fixed 
distance r between their central monomers. Field theoretical models for swollen polymer 
systems predict that the canonical partition function of a star polymer with k outgoing 
arms of length L is given by [8] : 

Zf ~ / L L fffc+fc<Tl (4) 

in the limit L — > oo. In Eq. (pH) \i is a non- universal parameter (the connectivity 
constant), while cr^'s are universal exponents associated to a vertex with k outgoing 
legs. In a star there are /c-ends (each contributing with an exponent <Ti) and one vertex 
with k outgoing legs (contributing with an exponent o^.), as shown in Fig. [1] (a). The 
value of the exponents is known from numerical simulations [H] and field theoretical 
computations using 4 — e expansions and resummation techniques [15] . 
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Figure 1. (a) An example of star polymer with k arms. There are k vertices of order 
1 and one of order k. (b) Nucleation of two polymers from the middle monomer. The 
contact leads to the formation of a star polymer with four outgoing legs, (c) Nucleation 
of two polymers from end monomers. The result after nucleation is only one polymer 
of double length with 2 vertices of order 1. 



Consider now two SAWs of length L, whose central monomers are separated by a 
vector r from each other. The partition function is given by: 

Z(r) ~ fi 2L L 4cTl g (5) 

where g(x) is a scaling function and v the SAW metric exponent. When the two polymers 
are sufficiently far from each other r 3> L u , or x — > oo one has g(x) — > 1, as the partition 
sum in Eq. (jSJ) should reduce to (^saw) 2 ; where ZgAw ~ fi L L 2ai is the partition function 
of a single self-avoiding walk. For r <^ L u , Eq. (jSJ) should reduce to that of a 4- arm star 
polymer Z$ 2 (see Fig. [T] (b)), which imposes the following scaling behavior in the limit 
x -t 

g(x) ~ x~ a4/u . (6) 

The numerical value for the exponent of three dimensional self-avoiding walks are 
o"4 ~ —0.48 pTTj . v ~ 0.588 [20], therefore —a^/u ~ 0.82. From the limiting behaviors 
for |r| — > oo and |r| — > a (a microscopic distance) one obtains the free energy 
difference between a configuration in which the two central monomers are far apart 



A lattice polymer study of DNA renaturation dynamics 



5 



and a configuration where the monomers are at a distance a as 

AF = k B Ta 4 (log L - - log a] (7) 



v 

One can repeat the same reasoning for contacts between end monomers (see 
Fig. [1] (c)). Let s be the separation between two end points of the polymers. The 
partition function is 

Z®~^L^h( l 4) (8) 



where again h(x) — > 1 for x — > oo, while the short distance limiting behavior is 

h{x) ~ x 2<T1 ^ (9) 

with o-i ~ 0.08, hence 2ai/V ~ 0.27. 

Similar scaling arguments were applied to the study of thermal denaturation of 
double stranded DNA 



3. Model 



We consider two self-avoiding walks of length L on a face centered cubic (fee) lattice, 
in which each site has 12 nearest neighbors displaced by (±1,±1,0), (±1,0, ±1), and 
(0,±1,±1). A configuration of the two walks is specified by the position of their 
monomers fi(i), ^(i), with i = 1, 2, . . . , L. The two walks are also mutually avoiding, 
except at homologous monomers which can overlap, i.e. ri(i) = r^C?) is only allowed if 
i = j. In the renaturated state the two walks overlap over their full length, fi(i) = ^(z) 
for all i. The version of this model on the simple cubic lattice was introduced by Causo 
et al. [6] for studying the thermodynamics of DNA denaturation, i.e. the unbinding of 
a double stranded polymer due to a temperature increase. 

In our simulations the starting configuration consists of two separate strands placed 
at random in a cubic box with periodic boundary conditions. The walks are first 
equilibrated while keeping their central monomer fixed. At time t = the polymers 
begin their diffusive motion until nucleation, i.e. the contact between homologous 
monomers, occurs. At that point the simulation ends and the nucleation time is 
registered. The polymer dynamics consists of local corner-flip or end-flip moves that 
do not violate self- avoidance. The probability that a given monomer moves is taken 
proportional to the number of nearest neighbors that can be reached on a fee lattice: 
we randomly attempt on average 11 Monte Carlo moves for each end-monomer and 3 
moves for each other monomer. 

We performed simulations for polymers of various lengths up to L = 193, within a 
cubic simulation box with side 100 (the lattice parameter being a = a/2). For the longest 
polymers considered in this work, the gyration radius is R g ~ 12, which is still much 
smaller than the box size, thus the simulations are performed in the diluted regime. We 
verified that the simulations reproduce the expected Rouse dynamics behavior D ~ 1/L 
[7J [33] for a self avoiding walk. This model neglects hydrodynamics effects, which 
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however are not expected to change the value of the exponents, at least in some regime 
(see Discussion). 

The Monte Carlo dynamics satisfies detailed balance, therefore it would be possible 
also to infer the renaturation rate from an initial bound state and waiting until unbinding 
takes place. Here we prefer to use the direct simulation approach, as it is closer to the 
experimental setup. 

4. Equilibrium properties 

We focus first on the equilibrium probability distributions of two strands to verify the 

scaling form of the partition functions derived in Section 3. We ignore nucleation events 

and in a single long run we sample the equilibrium configurations of the diffusing strands. 

Let us denote with r the distance between the middle monomers of the two strands and 

let P(r) be the probability distribution function for r. If r 3> R g , with R g the radius of 

gyration of the two polymers, then P(r) ~ n(r), where n(r) is the number of points at 

distance r from a given point in an fee lattice. Asymptotically in r we have n(r) ~ r 2 , 

equivalent to uniform probability in three dimensional space. 

The probability of finding the two central monomers at distances r < R g is 

suppressed due to mutual avoidance between the two interpenetrating strands. In order 

to get a good statistics on P(r) for r < R g we bias the dynamics so that configurations 

with small r are favored. This biased dynamics works as follows. A move changing 

the distance between the middle monomers from r to r' is accepted with probability 

min{l, n/n'}, where n and n' are the number of fee lattice points at a distance r and r' 

from the origin, respectively. The biased dynamics produces a probability distribution 

P(r) 

P*(r) which scales as P*(r) ~ — — in the limit of large r. This construction is 

particularly useful because we have P*(r) ~ g(rL~ u ), with g(x) introduced in Eqs. ([5])- 
(E]), that is, we can investigate directly the scaling function g(x). 

A plot of a P*(r) for two strands with L = 97 is shown in Fig. [2] (circles). 
For short distances the data reproduce the power-law behavior discussed in Eq. 
For intermediate distances P*(r) is constant, a signature that the two polymers are 
sufficiently far apart to be considered non-interacting (P(r) ~ r 2 ), whereas for longer 
distances the probability decays because of the finite size of the simulation box. 

In Fig. [2] (squares) we also show a similar plot for P*(s), the reweighted distance 
between the end monomers. For the calculation we used again biased dynamics, where 
the bias is on the distance between two end monomers and it is constructed with similar 
rules as for the central monomers. The results are consistent with the power-law limiting 
behavior of Eq. (jUJ) . 

5. Renaturation rates from simulations 

We now turn to the renaturation dynamics. As mentioned in the Introduction, the 
rate limiting step in DNA renaturation is the formation of a nucleus of a few base 
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Figure 2. Log-log plot of P*{r) and P*(s), the reweighted equilibrium distribution of 
distances between middle monomers (black circles) and end monomers (red squares) 
respectively. The polymers simulated consist both of 97 monomers. The decay at 
long distances is due to the finite size of the cubic box, while at short distances self- 
avoidance interactions play a central role. A fit of the data at short distances gives 
p*(r) ~ rO.saio.os^ and p*^ _ ^±0.05^ both in agreement wit h the scaling form 

discussed in Eq. (|6|) and Eq. ([9]). 



pairs, from where a rapid zippering follows. In the simulation we focus on the time 
needed to reach a nucleation event and the simulation is stopped when nucleation occurs. 
The average nucleation time (i nuc i.) is calculated from typically about 10 3 independent 
realizations, i.e. every time randomly generating the starting configuration of the two 
strands. This translates into an average elapsed real time of about two weeks on an 
Intel(R) Core(TM)2 Quad CPU Q9650 @ 3.00GHz processor. The nucleation rate is 
then obtained by: 

k= (tnucl.)" 1 - (10) 

It has been observed experimentally that the two strands must overcome a free 
energy barrier of about 30/csJ/mol [18] to form an active nucleus and that the nucleus 
corresponds to a double helix of about 15 bases [HUE]. Hence, the diffusion process may 
bring frequently two complementary monomers in contact, but many of these encounters 
are not "productive" as they do not lead to the formation of an active nucleus from 
which zippering of the helix would follow. We model this feature by introducing a 
local nucleation probability: a move that would bring two homologous monomers in 
contact [say r\(z) = r 2 (i)] is accepted with probability p n , otherwise it is rejected and 
the simulation continues. Hence, low values of p n correspond to high nucleation barriers 
and consequently long simulation times. To avoid confusion we will refer to the latter 
as a local barrier. This is not related to the free energy barrier AF originating from 
the mutual self avoidance of the renaturating strands discussed in the Introduction and 
given in Eq. (J7|). We will show below that the estimates of the scaling exponents depend 




Figure 3. Log-log plot of fc mon as a function of the polymer length L. k mon is obtained 
by means of Eq. (|10[) from an average of first contact times of middle monomers over 
different realizations with the polymers being initially placed at random in the box 
(size=100). As a result fc mon - L"0.58±o.oi for pn = i and ^ _ L -o.47±o.o4 for 

Pn = 0.1. 



on p n . 

Nucleation times corresponding to different values of p n can be obtained during 
the same run. Consider a set of nucleation probability values p n ^ < Pn < . . . < Pn\ 
each representing a different system. At an attempt of a nucleation step we generate a 
random number r uniformly distributed in the unit interval, and we consider as nucleated 
all systems q for which r < p$. Their nucleation time is added to the corresponding 
statistics, and the simulation continues in order to sample the nucleation times of the 
remaining systems, until nucleation for the smallest value pn ^ has been reached. 



5.1. Rate of middle monomers nucleation 

We consider first the nucleation between middle monomers. In this case the nucleation 
time is defined as the time for first contact between the two central monomers i = 
(L — l)/2 (for convenience we considered polymers with an odd number of monomers 
so that there is a unique middle monomer). In the course of the simulation two other 
analogous monomers may come in contact, say f\{k) = f 2 (k) with k ^ (L — l)/2, but 
the event is ignored and the simulation continues. 

In Fig. |3] we show our numerical results using nucleation probabilities p n — 1 and 
p n = 0.1 and for polymers of length up to L = 193. For each length we typically 
sample about 10 3 independent nucleation events. Both nucleation probabilities give a 
nucleation rate decaying as a power-law as a function of the length L. However the 
exponent depends on p n . Our estimates, which are based on linear fits for L > 30, 
yield k mon ~ £-o.58±o.oi f or ^ _ ^ (corresponding to the case of no local barrier) and 
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Figure 4. Log-log plot of fc 2 as a function of the polymer length L for different contact 
probabilities p n . The points are obtained by taking the inverse of the average time 
until first contact. Four series are shown corresponding to p n = 1, 0.1, 0.01, 0.001. 
The fit for Pn = 0.001, L > 30, gives k 2 ~ L°- 52±0,01 . 



k mon ~ L- - 4 ™ 4 for p n = 0.1. 

Kramers' theory (Eq. (j2])) [T7] predicts for A; mon a decay governed by an exponent 
<7 4 = —0.48. We note that our numerical data are in agreement with that for p n = 0.1. 
The results for p n — 1 are however not compatible with the Kramers' theory. 

5.2. Rate of polymer nucleation 

We consider next the nucleation rate k<i of the whole polymer. We define the nucleation 
time t nuc i as the first contact time between any of the L complementary monomers. 
Figure H] shows plots of k 2 as a function of L for different values of the nucleation 
probability (p n = 1,0.1,0.01,0.001). In absence of a local nucleation barrier (p n = 1) 
we find a non-monotonic behavior for It decreases for L < 50, while it increases 
for larger L. As the value of p n is lowered, k 2 becomes monotonic in L. From a fit 
of the data for L > 30 we obtain the estimates a = 0.34 ± 0.01 for p n = 0.01, while 
a = 0.52 ± 0.01 for p n = 0.001. Again, the exponent is consistent with that predicted 
by the Kramers' theory (Eq. (JH])) only for sufficiently low p n . 

5.3. Rate of end monomers nucleation 

We finally studied the renaturation dynamics for end monomers. In this case we define 
the nucleation time as the time needed to nucleate the monomers at one edge of the 
strand, at position i = 1. Figure [3] shows the data of fc en d as a function of L for two 
different values of the local nucleation probability. For p n = 1 we find as best fit to 
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Figure 5. Log-log plot of fc cn d as a function of the polymer length L. The points are 
obtained by taking the inverse of the average time until (on average) the 10 th contact 
between end monomers. The fit for L > 30 gives fc cn d ~ £-o.26±o.oi f or ^ _ ^ an( j 

fcend ~ L- 018±0 02 for Pn = 0.1. 



the data k cnd ~ L-°- 26±0 - 01 , while A; end ~ L-o.is±o.o2 for p n = 0.1. We also considered a 
smaller value of p n = 0.01 and found k cnd ~ £-0-i5±0.02 ^ata no i s hown). This suggests 
that a limiting value of the exponent has been reached. Kramers' theory would predict 
a scaling for the renaturation rate of the end monomers as k cnd ~ L~ 2<Tl = L~ 016 (see 
Eq. [3]). Again for sufficiently low p n (high local nucleation barriers) the simulations 
results are in agreement with Kramers' theory. 

6. Discussion 

In this paper we investigated a simple polymer lattice model of DNA renaturation. We 
focused on the scaling of the renaturation rate as a function of the polymer length. 
Three different rates were considered: the rates of renaturation for middle monomers 
{kmon), for end monomers (k en d) and for the whole strands (fo). We also introduced 
a local nucleation probability p n which is associated to an additional local free energy 
barrier. In DNA this barrier corresponds to the free energy one has to overcome to form 
an active nucleus of about 10-15 bases from where a rapid zippering starts. 

The simulations show that if the local barrier is sufficiently high (p n 1) the 
exponents governing the scaling behavior of k mon , k en d and &2 are all in agreement with 
Kramers' rate theory. An increase of the local barrier corresponds to an increase of 
the simulation time since the polymers attempt an increasing number of unsuccessful 
nucleation events. The process is slow, therefore the probability distribution of finding 
two monomers at a given distance from each other approaches its equilibrium value. 
In this regime Kramers' theory is expected to be valid. In addition, the whole strand 
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Figure 6. Schematic representation of the suggested experiment to measure the 
nucleation rate of end- monomers k en d- The strands have sequences of length I that are 
complementary at one end (labeled as 1 and 2). The rest of the strands is composed by 
a non-complementary L—l stretch (such as poly (A) stretches) so that nucleation cannot 
occur in this part. The end-renaturation dynamics can be detected by measuring the 
light emitted by a couple of quencher-fluorophore, as used in molecular beacons. 



renaturation rate k 2 scales with an exponent consistent with experimental data [T7] . 
which suggests that experimental conditions correspond to the high local barrier regime. 
Experiments so far focused on the rate k 2 , which is usually measured by monitoring the 
time dependence of the UV absorption spectrum. We argued, and verified in simulations 
in the high local barrier regime (p n <C 1), that the rate for end monomer nucleation 
scales as k cnd ~ L~ 018 . By appropriately engineering the sequences of the renaturating 
strands one could mesure /c cnc j in experiments. For this purpose it is sufficient to 
consider two sequences which are complementary in I end bases and non-complementary 
for the remaining L — I bases (see Fig. [6]). As hybridization can only occur at the 
edges, the dynamics is governed by k en d. The detection can be done by a quencher- 
fluorophore pair as used in molecular beacons [191 C] • The scaling of k en d is determined 
by increasing L while keeping I fixed. A modification of the sequence architecture, 
e.g. using complementary stretches in the middle of the sequence, would allow one to 
determine the renaturation rate for a single monomer k mon . 

We also performed simulations in the low local barrier regime (p n — > 1). In this case 
the exponents for k mon , k en d and k 2 deviate from Kramers' rate theory predictions. The 
renaturation is "fast" so that two complementary monomers in the approaching strands 
do not feel additional local barriers. In this regime it is natural to expect that at least one 
of the hypothesis behind Kramers' theory breaks down, namely the assumption that the 
two approaching strands can be considered as being in quasi-equilibrium conformations. 
It would be tempting to think that renaturation becomes a diffusion-limited reaction in 
this regime, but this is actually not true as there is still a free energy barrier (Eq. ([?])) 
generated by the self-avoidance interactions between the strands. A better quantitative 
understanding of this low local barrier regime remains challenging. We also expect that 
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the exponents in this fast renaturation regime would be influenced by hydrodynamics 
effects which are not built in in our simulations. This is an important point to bear in 
mind in view of possible comparison with experimental data. Hydrodynamics should 
not influence however the p n < 1 regime as renaturation is slow and occurs between 
strands in quasi-equilibrium regime. 

Another remark concerns the limit of validity of Eq. ([3]) . As the rate k<i grows with 
L (due to the increasing number of nucleation sites as L grows) we expect a deviation 
from Eq. ([3]) at some length L* beyond which renaturation becomes "too fast" in the 
sense discussed above: in this regime conformations of the approaching polymers are 
far from equilibrium. The scaling of &2 a s described by Kramers' theory should hold for 
L<L*. 

In conclusion, our simulations of renaturation show that the process has a very 
rich dynamical behavior. A key quantity appears to be the local free energy barrier 
for monomer-monomer nucleation. Different regimes appear in the limit of high barrier 
and of low barrier. The former regime is well-understood within the framework of 
Kramers' rate theory, while a better analytical understanding of the latter one remains 
challenging. 
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